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A  computational  study  is  performed  for  Mach  5  air  flow  over  a  cylinder  with  a  dielectric 
barrier  discharge  actuator  set  into  the  cylinder  surface.  The  actuator  is  pulsed  at  nanosec¬ 
ond  time  scales,  which  rapidly  adds  energy  to  the  flow,  thereby  creating  a  shock  wave 
that  travels  away  from  the  pulse  source.  As  the  shock  wave  travels  upstream,  it  interacts 
with  the  standing  bow-shock  and  momentarily  increases  the  bow-shock  standoff  distance. 

This  phenomenon  is  also  observed  in  phase-locked  schlieren  photography  captured  dur¬ 
ing  the  experiment.  The  focus  of  this  paper  is  to  reproduce  flow  phenomena  observed 
in  the  experiment  using  high-fidelity  computations  in  order  to  provide  additional  insight 
into  the  shock-shock  interaction,  the  effect  the  dielectric  barrier  discharge  pulse  has  on  the 
surface  properties  of  the  cylinder,  and  develop  a  reduced-order  phenomenological  model 
representative  of  the  nanosecond  pulse  discharge  system.  Experimental  and  high-fidelity 
modeling  studies  of  the  nanosecond  pulse  dielectric  barrier  discharge  plasma  actuators  are 
known  to  operate  with  relatively  low  temperatures.  This  work  explores  the  possibility 
that  the  induced  compression  wave  is  generated  by  rapid  thermalization  of  the  discharge 
which  results  in  a  local  temperature  rise  occurring  on  longer  time  scales.  Two-dimensional 
simulations  are  performed  and  provide  many  useful  details  about  the  discharge  event  while 
comparing  with  many  measurements  captured  by  the  experiment.  However,  the  simulations 
indicate  the  experiment  experiences  significant  three-dimensional  effects,  thus  requiring  a 
three-dimensional  simulation  of  the  entire  experiment  to  accurately  capture  the  complex 
cylinder/tunnel-sidewall  interaction  and  replicate  the  resultant  flow.  Three-dimensional 
results  of  the  discharge  event  reveal  the  discharge  pulse  produces  a  compression  wave  that 
interacts  with  the  standing  bow-shock  and  that  the  momentary  increase  in  the  bow-shock 
standoff  distance  is  not  due  to  the  interaction  the  three-dimensional  compression  wave  has 
with  the  cylinder /tunnel-sidewall  boundary  layer. 

Nomenclature 

equatorial  radii  of  an  ellipsoid,  [m] 
coefficient  of  pressure,  i/f~p°^2 
nondimensional  heating  coefficient,  y2  p”  u3 
Mach  number 
pressure,  [Pa] 
gas  constant,  [J/kg-K] 
heat  flux,  [W/m2] 

*  Research  Aerospace  Engineer,  High-Speed  Flow  Research  Group,  Computational  Aerophysics  Branch.  Member  AIAA. 

^Team  Lead,  High-Speed  Flow  Research  Group,  Computational  Aerophysics  Branch.  Associate  Fellow  AIAA. 

J Postdoctoral  Researcher,  Nonequilibrium  Thermodynamics  Laboratories,  Dept,  of  Mechanical  Engineering,  Member  AIAA. 

§  Professor,  Nonequilibrium  Thermodynamics  Laboratories,  Dept,  of  Mechanical  Engineering,  Associate  Fellow  AIAA. 


1  of  30 

This  material  is  declared  a  work  of  the  U.S.  Government  and  is  not  subject  to  copyright  protection  in  the  United  States. 

American  Institute  of  Aeronautics  and  Astronautics 


a,  6,  C 

C  p 

C  h 


Q 


Report  Documentation  Page 


Form  Approved 
OMB  No.  0704-0188 


Public  reporting  burden  for  the  collection  of  information  is  estimated  to  average  1  hour  per  response,  including  the  time  for  reviewing  instructions,  searching  existing  data  sources,  gathering  and 
maintaining  the  data  needed,  and  completing  and  reviewing  the  collection  of  information.  Send  comments  regarding  this  burden  estimate  or  any  other  aspect  of  this  collection  of  information, 
including  suggestions  for  reducing  this  burden,  to  Washington  Headquarters  Services,  Directorate  for  Information  Operations  and  Reports,  1215  Jefferson  Davis  Highway,  Suite  1204,  Arlington 
VA  22202-4302.  Respondents  should  be  aware  that  notwithstanding  any  other  provision  of  law,  no  person  shall  be  subject  to  a  penalty  for  failing  to  comply  with  a  collection  of  information  if  it 
does  not  display  a  currently  valid  OMB  control  number. 


1.  REPORT  DATE 

JAN  2012 


2.  REPORT  TYPE 


3.  DATES  COVERED 

00-00-2012  to  00-00-2012 


4.  TITLE  AND  SUBTITLE 

Computational  and  Experimental  Analysis  of  Mach  5  Air  Flow  over  a 
Cylinder  with  a  Nanosecond  Pulse  Discharge 

6.  AUTHOR(S) 


7.  PERFORMING  ORGANIZATION  NAME(S)  AND  ADDRESS(ES) 

Air  Force  Research  Faboratory,Wright-Patterson  AFB, OH, 45433-7512 

9.  SPONSORING/MONITORING  AGENCY  NAME(S)  AND  ADDRESS(ES) 


5a.  CONTRACT  NUMBER 


5b.  GRANT  NUMBER 

5c.  PROGRAM  ELEMENT  NUMBER 

5d.  PROJECT  NUMBER 


5e.  TASK  NUMBER 


5f.  WORK  UNIT  NUMBER 

8.  PERFORMING  ORGANIZATION 
REPORT  NUMBER 


10.  SPONSOR/MONITOR'S  ACRONYM(S) 

11.  SPONSOR/MONITOR'S  REPORT 
NUMBER(S) 


12.  DISTRIBUTION/AVAILABILITY  STATEMENT 

Approved  for  public  release;  distribution  unlimited 


13.  SUPPLEMENTARY  NOTES 


14.  ABSTRACT 


15.  SUBJECT  TERMS 

16.  SECURITY  CLASSIFICATION  OF: 


a.  REPORT 

unclassified 


b.  ABSTRACT 

unclassified 


c.  THIS  PAGE 

unclassified 


17.  LIMITATION  OF 

18.  NUMBER 

ABSTRACT 

OF  PAGES 

Same  as 

30 

Report  (SAR) 

19a.  NAME  OF 
RESPONSIBLE  PERSON 


Standard  Form  298  (Rev.  8-98) 

Prescribed  by  ANSI  Std  Z39-18 


Q 

T 

u,  v,  w 

U 

x,y,z 

a 

7 

£ 

e 

A 

p 

p 

Cfo 

T 

Subscript 

c 

cy 

9 

ve 

w 

0 

00 


=  total  energy  input  by  actuator,  [W] 

=  temperature,  [K] 

=  streamwise,  transverse,  and  spanwise  velocity  components,  [m/s] 

=  velocity  magnitude,  [m/s] 

=  streamwise,  transverse,  and  spanwise  coordinates 
=  accommodation  coefficient  for  partial  slip  wall  boundary  conditions 
=  specific  heat  ratio,  (1.4  for  air) 

=  emissivity 


=  angle  along  cylinder  surface,  tan  1 
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y-Vcy  \ 
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=  mean  free-path,  ^ [nr] 

=  kinetic  viscosity,  [kg/m-s] 

=  density,  [kg/nr3] 

=  Stefan-Boltznrann  constant,  5.6704  •  10~8  J/m2-s-K4 
=  duration  of  the  energy  deposition  pulse,  [s] 


=  center  of  ellipsoid 
=  center  of  the  cylinder 
=  gas 

=  vibrational-electron-electronic 
=  wall 
=  stagnation 
=  freestream 


I.  Introduction 

A  recent  experimental  study  of  a  nanosecond  pulse  Dielectric  Barrier  Discharge  (DBD)  in  a  Mach  5 
flow  demonstrated  the  feasibility  of  plasnra-based  supersonic  flow  control.1  A  bow-shock  perturbation  on  a 
microsecond  time  scale  was  detected  in  phase-locked  schlieren  visualization.  Compression  wave  generation 
due  to  a  rapid  localized  heating  from  the  DBD  propagates  upstream  from  the  cylinder  surface  and  interacts 
with  the  standing  bow-shock.  This  interaction  temporarily  increases  the  shock  standoff  distance.  This  series 
of  events  can  be  repeated  at  time  spacing  up  to  10  ps  (100  kHz).  Previous  demonstration  of  the  nanosecond 
pulse  DBD  includes  separated  flow  reattachment2  up  to  Mach  =  0.85,  characterization  of  compression  wave 
propagation  in  a  quiescent  air,3  and  visualization  of  large-scale,  spanwise  vortex  over  the  airfoil4  at  M  =  0.3. 

The  flow  control  mechanism  (rapid  heating)  in  these  experiments2, 4  appears  consistent  with  a  Localized 
Arc  Filament  Plasma  Actuator  (LAFPA).5-'  The  main  idea  of  this  approach  is  forcing  the  flow  with  a  high 
amplitude,  high  bandwidth  perturbation,  at  a  frequency  approaching  one  of  the  flow  instability  frequencies, 
thereby  triggering  subsequent  growth.  Previous  flow-control  studies  using  LAFPA  actuators  in  atmospheric 
pressure  jet  flows5-7  (M  =  0. 9-2.0),  demonstrated  significant  localized  heating  and  repetitive  shock- wave 
formation  by  the  plasma,  large-scale  coherent  structure  generation,  and  mixing  enhancement.  This  effect 
was  achieved  at  a  low  actuator  power  (~  10  W  per  actuator),  at  forcing  frequencies  near  the  jet  column 
instability  frequency  (preferred  mode).  This  low  power  budget  contrasts  with  previous  bow  shock  control 
studies,  typically  on  the  order  of  10  kW,  using  pulsed  DC  discharge,8  pulsed  microwave  discharge,9  and  laser 
optical  breakdown.10’ 11 

In  this  paper,  the  demonstrated  effect  will  be  reproduced  using  the  LeMANS  code,12,13  developed  at  Uni¬ 
versity  of  Michigan.14, 15  LeMANS  was  previously  used  for  hypersonic  flow-control  by  energy  deposition,12 
therefore  rapid  energy  coupling  by  the  nanosecond  pulse  DBD  is  modeled  in  the  same  phenomenological 
approach.  It  is  assumed  that  the  phenomenological  model  prediction  replicates  the  thermal  effect  of  the 
DBD.  Various  important  parameters,  such  as  temperature  distribution,  pressure  and  heat  coefficients  are 
not  available  in  the  experiment.  Thus,  this  modeling  study  is  meant  to  reveal  details  about  the  flow  per¬ 
turbation  mechanism,  and  afford  an  evaluation  of  its  validity  for  practical  applications,  such  as  flow  control 
at  hypersonic  inlet,  isolator,  and  engine  exhaust,  while  concurrent  work  by  Poggie  et  cil.16  focuses  on  high- 
fidelity  modeling  of  only  the  dielectric  barrier  discharge  in  order  to  identify  the  relevant  energy  transfer 
processes  and  develop  a  reduced-order  model  of  the  DBD. 

The  nonequilibrium  plasma  wind  tunnel  used  in  the  experiment  was  initially  manufactured  to  study 
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nonequilibrium  hypersonic  flows  and  develop  optical  diagnostics  that  can  be  portable  to  and  perform  mea¬ 
surements  at  national  ground  hypersonic  test  facilities.  Owing  to  this  purpose,  the  freestream  velocity 
and  temperature  have  been  measured.  Results  from  previous  measurements  of  the  freestream  velocity  and 
temperature  are  used  for  freestream  uncertainty  quantification. 

The  work  uses  LeMANS  to  compute  flow  around  a  5  mm  cylinder  and  a  6  mm  cylinder  using  two- 
dimensional  computations  to  compare  with  existing  experiments.  While  the  5  mm  cylinder  simulation 
results  match  the  experiment,  the  6  mm  cylinder  simulation  results  over-predict  the  bow-shock  location. 
The  nanosecond  DBD  experiment  used  the  6  mm  cylinder,  so  a  parametric  study  exploring  the  uncertain¬ 
ties  associated  with  the  freestream,  wall  boundary,  and  thermo-chemical  nonequilibrium  flow  conditions  is 
conducted  to  quantify  their  effects  on  the  standing  bow-shock  location.  The  study  revealed  that  none  of  the 
parameters  explained  the  discrepancy  in  bow-shock  location.  Nonetheless,  the  two-dimensional  simulations 
allow  the  large  design-space  to  be  explored  with  minimal  computational  cost,  including  results  from  the  phe¬ 
nomenological  energy  deposition  modeling  of  the  discharge  event.  The  energy  deposition  modeling  results 
are  able  to  replicate  the  DBD  discharge  induced  compression  wave  speed  and  resultant  perturbed  bow-shock 
shape  and  qualitatively  match  phase-locked  schlieren  images. 

In  order  to  explain  the  bow-shock  standoff  discrepancy,  a  three-dimensional  computation  of  the  full 
tunnel  is  computed  and  found  to  agree  well  with  all  available  experimental  measurements.  The  results  show 
that  cylinder/tunnel  sidewall  interaction  produces  a  significant  pressure  drop  along  the  cylinder  span,  which 
draws  the  standing  bow-shock  closer  to  the  cylinder.  In  addition,  a  three-dimensional  simulation  with  the 
phenomenological  energy  deposition  model  representing  the  discharge  event  is  able  to  reproduce  many  of  the 
features  observed  in  the  phase- locked  schlieren  images  of  the  experiment. 

II.  Experimental  Facilities 

A  schematic  of  a  small-scale  Mach  5  nonequilibrium  wind  tunnel  is  shown  in  Fig.  1  (also  used  in  Ref. 
17, 18).  The  wind  tunnel  was  operated  using  dry  air  supplied  from  high-pressure  cylinders,  at  plenum 
pressures  of  po  =  370  Torr  (0.5  atm),  and  the  mass  flow  rate  of  7  g/s.  The  steady  state  run  time  at  the 
constant  static  pressure  in  the  supersonic  test  section  is  up  to  10  seconds.  The  flow  expands  through  an 
aerodynamically  contoured  Mach  5  nozzle,  with  the  throat  height  of  1.6  mm.  Top  and  bottom  walls  of  the 
supersonic  test  section  after  the  nozzle  exit  each  diverge  at  a  1.5°  angle  to  provide  boundary  layer  relief. 
Rectangular  optical  access  windows,  made  of  UV-grade  fused  silica,  are  flush  mounted  in  all  four  walls  in 
the  supersonic  test  section,  providing  ample  optical  access  for  schlieren  photography,  emission  spectroscopy, 
Planar  Laser  Induced  Fluorescence  (PLIF)  visualization  /  thermometry,17, 18  and  NO2  Molecular  Tagging 
Velocimetry  (MTV)19  diagnostics,  as  well  as  currently  on-going  vibrational  temperature  measurement  by 
picosecond  Coherent  Anti-Stokes  Raman  Scattering  (CARS)  system.20 


Diffuser 


Figure  1.  Schematic  of  the  Mach  5  wind  tunnel  with  a  cylinder  model  installed  in  the  test  section. 


The  dimensions  of  a  quartz  cylinder  model  used  to  generate  a  bow-shock  and  also  used  as  a  DBD  actuator 
is  4  cm  long,  6  mm  outer  diameter,  and  4  mm  inner  diameter.  The  model  is  located  14.5  cm  downstream 
of  the  throat  (3.5  cm  downstream  of  the  end  of  nozzle  contour),  where  the  flow  cross  sectional  area  is  4 
cm  x  4.6  cm.  The  ends  of  the  model  are  embedded  in  the  optical  access  windows  in  the  side  walls  of  the 
test  section.  The  baseline  shock  standoff  distance  (without  nanosecond  pulse  discharge),  measured  in  the 
schlieren  photography  is  1.2  mm,  with  the  spanwise  length  of  about  1  cm,  which  is  25%  of  the  test  section 
width.  The  spanwise  extent  of  the  bow-shock  is  limited  by  the  boundary  layer  growth  on  the  sidewalls  of 
the  test  section.17 
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The  test  section  static  pressure,  p^  =  1.2  Torr  (160  Pa),  was  measured  using  a  wall  pressure  tap  in  the 
side  wall  at  the  end  of  the  nozzle,  and  4  cm  upstream  of  the  cylinder  model.  The  pressure  is  assumed  to  be 
constant  through  the  sidewall  boundary,  a  valid  assumption  for  a  laminar  boundary  layer  (no  experimental 
measurements  taken  thus  far  have  indicated  the  wall  boundary  layer  is  turbulent).  The  flow  Mach  number 
inferred  from  the  plenum  pressure  and  freestream  static  pressure  is  =  4.6,  while  the  total  pressure 
measured  at  the  flow  stagnation  point  (downstream  of  the  bow-shock)  was  p  =  36  Torr,  which  corresponds 
to  a  freestream  Mach  number  of  Moo  =  4.8,  using  one-dimensional  normal  shock  relations  and  assuming  the 
static  pressure  is  constant  across  the  side-wall  boundary  layer. 

Figure  2  shows  a  cartoon  of  the  electrode  configuration  for  the  bow-shock  perturbation  by  the  nanosecond 
pulse  DBD.  One  actuator  electrode  is  composed  of  a  1  cm  long,  3  mm  diameter  copper  tube  with  a  tube  wall 
thickness  of  0.35  mm,  and  is  immersed  inside  the  quartz  tube.  Since  the  tube’s  inside  diameter  is  4  mm,  the 
actuator  electrode  is  positioned  to  contact  the  quartz  surface  near  the  upstream  side,  as  seen  in  Fig.  9.  The 
other  electrode  consists  of  a  strip  of  adhesive  copper  tape  1.5  mm  wide  and  12  mm  long,  attached  to  the 
outside  surface  of  the  quartz  tube  model.  The  two  electrodes  overlap  over  a  spanwise  distance  of  10  mm, 
centered  behind  the  core  flow  bow-shock,  with  the  ends  of  the  tape  covered  by  non-conducting  Kapton® 
tape.  Output  pulse  voltage  and  current  were  measured  during  each  run.  Measured  peak  voltage  and  current 
are  27  kV  and  70  A,  respectively,  coupling  between  4-7  mJ/pulse. 


Copper 

tape 


Positive 
polarity  pulse 


Negative 

polarity 

pulse 


Figure  2.  Diagram  of  the  cylinder  model  with  a  nanosecond  pulse  surface  DBD  plasma  actuator. 


Time  evolution  of  the  shock  generated  by  the  NanoSecond  Dielectric  Barrier  Discharge  (NSDBD),  as 
well  as  its  interaction  with  the  bow  shock  was  recorded  using  a  phase-locked  schlieren  system.3  It  was 
shown  that  a  discharge  pulse  generates  a  compression  wave  that  propagates  upstream  and  locally  ‘pushes’ 
the  bow-shock  away  from  the  cylinder.  This  perturbed  region  bends  away  from  the  flow  stagnation  line, 
convects  downstream,  and  eventually  returns  the  shock  envelop  to  the  baseline  shock  shape  about  20  p.s 
later.  Image  sets  for  this  microsecond-scale  shock-shock  interaction  were  taken  both  for  a  ‘single  pulse’  mode 
(pulse  repetition  rate  of  200  Hz),  and  ‘double  pulse’  mode  (two  pulses  separated  by  a  10  microsecond  delay, 
which  corresponds  to  a  pulse  repetition  rate  of  100  kHz).  This  study  will  primarily  use  on  the  ‘single  pulse’ 
mode,  so  each  discharge  pulse  acts  on  the  baseline  bow-shock.  Additional  details  about  the  NSDBD  process 
are  available  in  Ref.  1. 

A.  Freestream  flow  parameters 

The  facility  is  able  to  achieve  Mach  5  flows  using  a  blow-down  wind  tunnel.  While  it  is  not  possible  to 
exactly  characterize  the  test  section  for  each  experiment  run,  through  various  proven  approaches,  the  inputs 
necessary  for  a  computational  simulation  can  be  determined.  As  shown  in  previous  work,17  the  test  section 
has  an  inviscid  core  that  composes  about  25%  of  the  cross-section  area.  Since  the  cylinder  test  model  is 
essentially  two-dimensional  and  only  composes  a  fraction  of  the  core  flow,  it  is  reasonable  to  assume  that  the 
freestream  flow  is  uniform  and  the  resultant  flow-field  will  be  two-dimensional.  This  assumption  is  justified 
by  looking  at  Fig.  3,  a  top-down  view  of  phase-locked  schlieren  image  at  po  =  370  Torr  dry  air  (baseline 
condition),  with  the  standoff  distance  AS  =  1.2  mm. 

Freestream  velocity  was  measured  Uqq  =  719  ±  6  m/s  in  the  NO2  MTV.19  In  this  measurement,  a  small 
amount  of  NO2  (~  1%  mole  fraction),  was  seeded  in  the  main  nitrogen  flow  at  the  plenum  pressure,  p0  =  370 
Torr.  A  pulse-burst  laser  system21  generated  two  simultaneous  outputs:  one  for  the  photo-dissociation  of 
NO2  into  NO  and  O  at  355  nm  to  tag  a  line  in  the  flow,  the  other  tuned  at  226  nm  for  interrogation  of  the 
line  progression  by  NO  PLIF  imaging.  The  NO  fluorescence  intensity  distributions  were  fitted  by  a  Gaussian 
curve,  and  it  was  found  that  the  average  velocity  is  719  m/s  with  absolute  statistical  error  of  a  =  5.8  m/s. 
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Figure  3.  Phase-locked  schlieren  image  (top-down  view)  in  dry  air  at  po  =  370  Torr  at  baseline  condition 
(without  discharge).  Flow  direction  is  top  to  bottom  in  the  figure. 


The  NO  PLIF  thermometry  was  used  for  temperature  measurement  behind  the  Mach  5  bow-shock  in 
Ref.  18.  A  sheet  of  laser  was  tuned  to  pump  a  rotational  transition  of  J  =  5.5  or  J  =  16.5  on  the  NO(X,v'  = 
0  — >A,v"  =  0)  band.  The  pair  of  these  lines  was  selected  because  the  ratio  of  absorption  strength  is  linear 
at  around  To  =  300  K.  In  one  of  preliminary  experiments,  sheet  laser  wavelengths  were  tuned  to  excite  the 
J  =  5.5  and  J  =  9.5  rotational  states  in  order  to  probe  lower  temperatures  in  the  freestream.  Only  a  few  of 
experiments  have  been  carried  out  thus  far  with  some  uncertainty  in  the  results.  But,  the  results  tend  to  fall 
within  a  temperature  range  of  about  T^  =  50—60  K,  which  appears  consistent  with  a  freestream  temperature 
obtained  by  considering  the  freestream  pressure  ratio,  T^  =  56  K.  Freestream  temperatures  calculated  from 
the  isentropic  relation,  using  the  freestream  velocity  measured  in  the  NO2  MTV  measurement,  are  Too  =  78 
K,  60  K,  and  50  K  for  Mach  numbers  of  M  =  4,  4.55,  and  5,  respectively.  The  freestream  temperature  for 
M  =  4.55  is  also  close  to  this  isentropic  flow  temperature. 

As  previously  mentioned,  the  test  section  static  pres¬ 
sure,  Poo  =  1.2  Torr  (160  Pa),  was  measured  using  a  wall 
pressure  tap  in  the  side  wall  at  the  end  of  the  nozzle,  and 
4  cm  upstream  of  the  cylinder  model.  The  precision  on 
the  pressure  gauge  is  0.1  Torr,  so  its  uncertainty  is  ±0.05 
Torr  (±6.67  Pa).  The  freestream  density  is  inferred  using 
the  ideal  gas  relation  (p  =  pRT).  For  this  work,  the 
freestream  dry  air  is  composed  of  78%  nitrogen  (N2),  and 
22%  oxygen  (02),  by  density.  Table  1  lists  the  nominal 
freestream  conditions  and  uncertainty  bounds.  The  effect 
of  the  uncertainties  will  be  addressed  in  a  later  section. 


Table  1.  Freestream  conditions  and  uncertainty 
bounds  for  Mach  5  air  flow  around  a  cylinder. 


Parameter 

Value 

Uoo,  [m/s] 

719  ±6 

Too,  [K] 

56  ±5 

Poo,  [kg/rn3] 

0.009924  ±  0.0013 

Poo,  [Pa] 

160  ±6.7 

Moo 

4.76  ±0.25 

III.  Numerical  Method 

Flow-field  results  are  obtained  using  Computational  Fluid  Dynamics  (CFD)  to  solve  the  Navier-Stokes 
equations.  The  CFD  computations  are  executed  using  the  Michigan  Aerothermodynamic  Navier-Stokes  (Le- 
MANS),  code  developed  at  the  University  of  Michigan.14, 15  LeMANS  is  a  general  2D/axisynnnetric/3D, 
parallel,  unstructured  finite-volume  CFD  code  and  has  been  used  previously  in  numerous  studies  of  hy¬ 
personic  flows. 12-15  LeMANS  may  be  employed  with  any  of  three  thermodynamic  models:  perfect  gas, 
equilibrium,  and  nonequilibrium  thermochemistry.  LeMANS  employs  a  two-temperature  model  to  account 
for  thermal  nonequilibrium  and  a  standard  finite-rate  chemistry  model  for  nonequilibrium  chemistry.  The 
two  temperature  model  assumes  that  a  single  temperature,  T,  accounts  for  the  translational  and  rotational 
energy  modes  of  all  species  while  the  vibrational  and  electronic  energy  modes  are  accounted  for  by  a  sep¬ 
arate  temperature,  Tve.  The  simulations  are  performed  using  second-order  accurate  spatial  discretization 
and  carry  double  precision  arithmetic  throughout. 

LeMANS  is  primarily  used  for  steady-state  simulations,  but  is  capable  of  computing  time  accurate  sce¬ 
narios  with  first-order  temporal  accuracy.  However,  numerical  error  associated  the  low  temporal  accuracy  is 
minimized  by  enforcing  a  sub-nanosecond  time  step  (At  <  1  ns). 

The  nanosecond  DBD  discharge  used  in  the  experiment  is  effectively  a  thermal  actuator.  As  such,  a 
phenomenological  model  of  dissipative  heating  is  used  to  represent  it.  This  model  is  accounted  for  in  the 
Navier-Stokes  equations  by  the  addition  of  a  source  term,  S,  to  the  right  side  of  total  energy  equation, 


5  of  30 


American  Institute  of  Aeronautics  and  Astronautics 


without  any  energy  being  directly  deposited  into  the  vibrational-electron-electric  energy  equation  when  the 
simulation  is  performed  assuming  thermodynamic  nonequilibrium.  Deposition  of  all  the  energy  into  the 
translational  mode  is  a  strong  assumption,  but  is  adequate  for  the  purpose  of  this  study  since  it  is  assumed 
the  compression  wave  generated  by  the  DBD  is  due  to  a  rapid  transfer  of  energy  into  the  translation  energy 
mode,  an  observation  seen  in  previous  work  by  Popov.22 

The  shape  and  location  of  the  actuator  are  modeled  with  contours  of  constant  S  having  an  ellipsoidal 
shape.  This  approach  has  been  used  successfully  in  previous  numerical  investigations.12  The  strength  (total 
energy),  deposited  into  the  flow  uses  exponential  decay  from  the  centroid  of  the  energy  deposition  pattern, 
which  for  two-dimensional  simulations  is: 

«=  St,  «p(-(§)2-(|)2) 

o  )  X  '  (1) 

X  =  [x  —  Xc) 

y  =  {y-  yc) 

where  variables  a  and  b  are  the  equatorial  radii  (along  the  x  and  y  axes).  Coordinates  (xc,  yc )  represent  the 
centroid  of  the  ellipsoid.  Note  that  Q  represents  the  total  power  deposited  in  the  flow  and  f Sdxdy  =  Q. 

A.  Grid  Independence 

Figure  4  shows  a  structured  mesh  for  a  two-dimensional  computational  domain  around  a  6  mm  diameter 
cylinder.  Flow  direction  is  left  to  right,  and  the  stagnation  point  is  located  at  x  =  0.003  m.  Only  the  first  90° 
of  the  cylinder  surface  is  computed  in  this  section  to  minimize  computational  cost  and  because  measurements 
behind  the  cylinder,  such  as  the  wake  or  shear  layer  structure,  were  not  captured  in  the  experiment.  Note  full 
cylinder  computations  are  discussed  in  subsequent  sections.  The  meshes  were  developed  such  that  clustering 
in  the  radial  direction  occurs  at  both  the  cylinder  surface  and  at  the  bow-shock,  while  the  mesh  in  the 
azimuthal  direction  is  distributed  to  provide  orthogonality  of  the  cells  at  the  bow-shock. 

0.009r 


X  [m] 


Figure  4.  Shock- fitted  grid  for  Mach  5  air  flow  around  a  cylinder.  (100  X  100) 

Three  doubly  refined  grids:  50  x  50  (coarse),  100  x  100  (medium),  and  200  x  200  (fine),  were  used  for 
the  grid  independence  study.  Input  parameters  were  taken  from  the  nominal  conditions  listed  in  Table 
1;  Uao  =  719  m/s,  T^  =  56  K,  and  =  0.009924  kg/m3.  Since  the  input  parameters  correspond  to 
the  baseline  scenario  (i.e. ,  steady-state),  the  solution  is  advanced  until  the  root-mean-square  residual  error 
approaches  machine-precision  and  remains  unchanged  for  subsequent  time  steps,  as  seen  in  Fig.  5.  As  the 
mesh  number  increases,  the  number  of  iterations  required  for  the  convergence  also  increases  due  to  the 
reduction  in  the  time  step  required  to  maintain  stability,  which  is  calculated  based  on  the  smallest  cells 
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size  used  in  the  computational  domain.  In  addition,  a  maximum  time  step  of  At  =  1  ps  is  enforced  for  all 
steady-state  simulations  performed  in  this  paper. 


Figure  5. 


Root-mean-square  residual  error  for  the  nominal  baseline  case  with  the  medium  grid. 


Figure  6  compares  the  temperature  contours  from  the  fine  grid  (top)  and  medium  grid  (bottom).  The 
distributions  are  nearly  identical,  except  at  the  bow-shock,  where  the  shock  appears  more  smeared  for  the 
medium  grid.  While  the  shock  thickness  decreases  with  increasing  mesh  number,  the  standoff  distance 
(defined  as  the  distance  between  the  stagnation  point  and  the  location  of  maximum  density  gradient), 
remains  the  same  (AS  =  1.55  mm),  for  all  three  grids. 


Figure  6.  Temperature  contours  for  Mach  5  air  flow  over  a  cylinder  with  fine  grid  (top)  and  medium  grid 
(bottom). 

Figure  7  plots  pressure  coefficient,  cp,  and  nondimensional  heating  coefficient,  c h,  along  the  cylinder 
surface  for  all  three  grids.  Coordinates  along  the  cylinder  are  converted  into  the  degree  angle,  9,  with 
0  =  0°  located  at  the  stagnation  point.  The  pressure  coefficient  is  nearly  identical  for  all  three  grids, 
while  the  nondimensional  heating  coefficient  exhibits  an  offset  for  the  coarse  grid.  Both  the  pressure  and 
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nondimensional  heating  coefficients  achieved  on  the  fine  grid  replicate  those  from  the  medium  grid.  Therefore, 
the  medium  grid  (lOOx  100),  is  considered  grid-independent  and  will  be  used  in  the  remaining  two-dimensional 
baseline  simulations,  unless  otherwise  noted. 


Figure  7. 


Pressure  and  nondimensional  heat  coefficient  along  the  cylinder  surface  for  various  grids. 


The  temperature  distribution  is  compared  with  a  schlieren  image  of  a  bow-shock  ahead  of  6  mm  diameter 
cylinder  in  Fig.  8(a).  The  experiment’s  bow-shock  standoff  distance  of  AS  =  1.2  mm  is  20%  smaller  than 
the  computed  standoff  distance  of  AS  =  1.55  mm.  Also  displayed  in  the  figure  is  an  empirical  relation  for 
the  shock  envelop  developed  by  Billig  et  al.  23,24  which  was  developed  by  comparing  schlieren  images  from 
many  shock  tube  experiments  and  has  dependencies  on  Mach  number  and  cylinder  radius.  Billig’s  empirical 
relationship  to  determine  the  bow-shock  standoff  location  on  the  stagnation  line  is  specified  later  in  Eq.  (3) . 
For  a  Mach  number  of  M  =  4.8  and  the  radius  of  3  mm,  Billig’s  empirical  standoff  distance  AS  =  1.42  mm  is 
much  closer  to  the  LeMANS  computed  standoff  distance  and  shows  a  good  agreement  with  the  shock  profile 
as  well. 


(a)  6  mm 


(b)  5  mm 


Figure  8.  Contour  lines  of  temperature  for  Mach  5  air  flow  over  a  cylinder  for  two  cylinder  diameters.  The 
figures  include  schlieren  images  from  the  experiment  and  Billig’s  empirical  formula  of  the  shock  envelope.24 
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The  discrepancy  in  the  bow-shock  standoff  location  between  the  computation  and  the  experiment  sug¬ 
gests  there  may  be  inappropriate  assumptions  made  in  the  simulation  or  there  are  three-dimensional  effects 
in  the  experiment  that  cannot  be  replicated  by  the  two-dimensional  computation.  To  help  assess  the  validity 
of  the  two-dimensional  computational  results,  a  computation  was  performed  for  a  5  mm  diameter  cylinder 
scenario,  which  was  also  performed  experimentally  (although  the  nanosecond  DBD  was  not  employed  for 
this  experimental  setup).  As  seen  in  Fig.  8(b),  the  5  mm  cylinder  experimental  results  are  in  much  better 
agreement  with  the  computations  and  Billig’s  empirical  relationship.  These  results  suggest  that  the  differ¬ 
ences  in  standoff  distance  observed  in  the  6  mm  diameter  cylinder  scenario  are  mostly  likely  attributed  to 
a  three-dimensional  effect  found  in  the  experiment  or  greater  uncertainty  in  freestream  and  wall  boundary 
conditions.  The  following  subsections  explore  the  effects  associated  with  tunnel  and  boundary  conditions 
using  the  two-dimensional  simulations  since  exploration  of  the  full  three-dimensional  simulation  requires  a 
significant  increase  in  computational  cost. 

B.  Wall  Boundary  Conditions 

One  uncertainty  about  the  experiment  is  the  temperature  of  the  cylinder  surface.  Since  the  cylinder 
is  initially  at  room  temperature,  one  assumption  is  that  the  cylinder  surface  is  at  a  constant  isothermal 
condition.  However,  given  the  experiment’s  long  run  time  (10  seconds),  the  cylinder  should  cool  to  adiabatic 
conditions  (i.e. ,  no  heat  transfer  to  the  surface).  Another  common  approach  in  CFD  is  to  assume  the  surface 
is  in  radiative  equilibrium.  This  assumption  is  modeled  by  using  the  Stefan-Boltzmann  Law: 

qw=ea0  T*  (2) 

where  qw  is  the  heat  transfer  to  the  wall,  e  is  the  emissivity  of  the  body,  and  cto  is  the  Stefan-Boltzmann 
constant.  Assuming  the  body  (i.e.,  cylinder  surface),  can  be  treated  as  a  ‘black  body’,  its  emissivity  is  unity 
(e  =  1).  However,  Fig.  9  shows  the  cylinder  is  composed  of  quartz  (a  dielectric),  which  has  an  emissivity25 
of  s  =  0.93,  except  for  the  exposed  electrode,  which  has  an  emissivity26  of  e  =  0.02.  Preliminary  simulations 
accounting  for  the  actual  thickness  of  the  exposed  electrode  (0.1  mm),  extruded  from  the  quartz  surface  did 
not  significantly  alter  the  resultant  flow-field,  so  the  remaining  simulations  assume  the  cylinder  surface  has 
a  continuous  radius  of  3  mm  (i.e.,  the  exposed  electrode  is  incorporated  into  the  cylinder  surface).  Table  2 
tabulates  the  four  scenarios  considered  for  the  wall  boundary  condition  study.  Each  steady-state  scenario 
uses  the  shock- fitted  grid- independent  mesh  and  nominal  freestream  conditions  listed  in  Table  1,  without 
the  nanosecond  discharge  event  (i.e.,  baseline  scenario). 

Quartz  tube 

freesteam 

*  |  negative 

'w  *  pulse 

- ► 

positive 

pulse 


Figure  9.  Side-view  cartoon  of  the  experimental  setup  of  the  cylinder. 


Table  2.  Wall  boundary  condition  scenarios  investigated  for  Mach  5  air  flow  over  a  cylinder. 


Scenario 

Details 

Adiabatic 

o 

II 

Isothermal 

T™  =  300  K 

Fully  Radiative 

£  =  1 

Mixed  Emissivity 

^quartz  =  0.93 

^•copper  —  0.02 
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Figures  10(a)  and  10(b)  show  the  surface  temperature  distribution  and  the  nondimensional  heating 
coefficient,  along  the  surface  for  each  of  the  scenarios.  As  seen  in  Fig.  10(a),  the  adiabatic  solution  is  very 
similar  to  the  radiative  boundary  condition  scenarios,  whereas  the  nondimensional  heating  coefficient  is  very 
similar  for  the  isothermal  and  radiative  boundary  conditions. 


(a)  Surface  temperature. 


(b)  Nondimensional  heating  coefficient. 


Figure  10.  Distributions  along  the  surface  of  a  cylinder  in  Mach  5  air  with  various  wall  boundary  conditions. 


Table  3  tabulates  the  shock  standoff  distance,  which  is  de¬ 
termined  by  identifying  the  peak  density  gradient  in  stream- 
wise  direction  (i.e.,  along  the  x— axis).  While  the  surface  pro¬ 
files  seen  in  Figs.  10(a)  and  10(b)  exhibit  differences,  the  shock 
standoff  distance  is  not  influenced  by  the  wall  boundary  con¬ 
dition  selected.  This  occurs  because  the  boundary  layer  that 
develops  on  the  cylinder  surface  is  quite  small,  and  since  the 
adiabatic  wall  boundary  condition  only  significantly  influences 
the  boundary  layer  thickness,  the  shock  standoff  distance  is  not 
influenced  by  the  boundary  condition  enforced.  The  effect  of 
wall  boundary  conditions  is  expected  to  be  more  pronounced  at  higher  temperatures,  such  as  during  the 
nanosecond  discharge  event,  where  the  gas  temperature  near  the  stagnation  point  reaches  2000  K  on  a 
microsecond  time  scale. 

C.  Thermo-Chemical  Nonequilibrium 

The  freestream  conditions  are  quite  cold  (see  Table  1),.  As  such,  it  is  unlikely  that  a  significant  amount 
of  thermal-chemical  nonequilibrium  is  present  in  the  flow.  For  completeness,  with  regards  to  chemical 
nonequilibrium,  simulations  were  run  using  a  one-species  perfect  gas,  and  a  five-species  gas  (N2,  O2,  NO, 
N,  and  O),  using  Park’s  1990  data  sets.  Since  the  latter  of  the  two  scenarios  requires  LeMANS  to  carry  5 
conservation  equations  for  the  mass,  the  computational  resources  required  to  obtain  a  solution  also  increases. 
As  anticipated,  the  simulation  results  show  no  appreciable  differences  between  the  one-species  perfect  gas 
and  the  five-species  finite-rate  chemistry  simulations.  For  completeness,  two  simulations  are  also  conducted 
allowing  for  thermal  nonequilibrium  of  the  vibrational  temperature.  The  Landau-Teller  model2'  is  used 
to  account  for  energy  exchange  between  the  vibrational-electronic  and  the  translational-rotational  energy 
modes. 

Typically,  rates  of  vibration-translation  energy  transfer  in  a  nitrogen-oxygen  system  are  very  low,28  and 
therefore  vibrational  temperature  at  the  plenum  T„e=300  K  is  assumed  frozen  in  the  freestream.  For  this 
thermo-chemical  nonequilibrium  case,  translational  and  vibrational  temperature  distributions  are  compared 
in  Fig.  11(a).  For  comparison,  Fig.  11(b)  displays  translational  and  vibrational  temperature  distributions 


Table  3.  Shock  standoff  distance  for  various 
wall  boundary  scenarios. 


Scenario 

Shock  Standoff 

Adiabatic 

1.55  mm 

Isothermal 

1.55  mm 

Fully  Radiative 

1.55  mm 

Mixed  Emissivity 

1.55  mm 
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at  the  freestream  vibrational  temperature,  T„e=240  K,  which  corresponds  to  a  scenario  where  thermal 
equilibrium  exists  up  to  the  nozzle  throat,  at  which  point  it  freezes  and  remains  frozen  downstream  into  the 
test  section. 
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-0.003 
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0.006 


0.003 
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-0.003 


-0.006 
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(a)  Tllei00  =  300  K 


(b)  Tve, oc  =  240  I< 


Figure  11.  Translational  and  vibrational  temperature  contours  for  Mach  5  air  flow  over  a  6mm  cylinder. 


As  seen  in  the  figures,  the  vibrational  temperatures  are  essentially  frozen  across  the  shock.  Also,  no 
noticeable  differences  in  the  standoff  distance  or  wall  properties  are  observed  in  either  scenario.  Therefore, 
the  thermal  nonequilibrium  effects  are  ignored  for  the  remainder  of  the  baseline  calculations.  However, 
finite-rate  chemistry  is  included  in  the  energy  deposition  scenarios  due  to  the  rapid  rise  in  local  temperature 
during  the  discharge  event. 


D.  Freestream  Uncertainty  Quantification 

The  freestream  conditions  observed  in  the  wind  tunnel  have  uncertainty  associated  with  them.  As  such, 
it  is  important  to  quantify  the  effect  these  uncertainties  have  of  flow  properties  of  interest,  such  as  standoff 
distance,  as  they  may  explain  the  discrepancy  observed  for  the  6  mm  diameter  scenario.  LeMANS  requires 
the  freestream  density,  temperature,  and  velocity  as  input  parameters.  These  three  variables  determine 
the  possible  scenarios  that  need  to  be  simulated  to  quantify  the  uncertainty  in  baseline  shock  standoff 
distance  associated  with  the  nominal  freestream  conditions.  Table  4  lists  the  various  runs  with  the  parameter 
combinations.  The  simulations  were  performed  assuming  thermal-chemical  equilibrium  and  the  isothermal 
wall  boundary  condition,  T^,  =  300  K. 

Results  of  standoff  distance  are  summarized  in  Fig.  12,  which  plots  standoff  distance  versus  Mach  number. 
The  points  are  denoted  by  color  to  indicate  the  Reynolds  number  (based  on  cylinder  radius).  As  seen  in  the 
figure,  the  standoff  distance  decreases  with  increasing  Reynolds  number.  The  general  trend  of  these  points 
appears  to  be  consistent  with  Billig’s  empirical  formula  for  cylinder  flow.24 


AS 

r 


0.386  exp 


(3) 
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Table  4.  Input  parameters  for  uncertainty  quantification  of  freestream  conditions  for  Mach  5  air  flow  over  a 
cylinder. 


Run 

Moo,  [m/s] 

Too,  [K] 

Poe  [kg/m3] 

R e^/L  [in'1] 

Mach 

1 

713 

50 

0.008876 

1.97-  106 

5.02 

2 

713 

50 

0.011578 

2.57-  106 

5.02 

3 

713 

60 

0.008876 

1.59-  106 

4.59 

4 

713 

60 

0.011578 

2.07-  106 

4.59 

5 

725 

50 

0.008876 

2.0  •  106 

5.11 

6 

725 

50 

0.011578 

2.61  •  106 

5.11 

7 

725 

60 

0.008876 

1.62  •  106 

4.66 

8 

725 

60 

0.011578 

2.11  •  106 

4.66 

where  AS  is  the  shock  standoff  distance,  r  is  the  cylinder  radius,  and  M  is  the  upstream  Mach  number. 
The  similarity  between  Eq.  (3)  and  the  solutions  in  Fig.  12  exists  because  expansion  of  Eq.  (3)  generates 
the  quadratic  dependence  of  the  standoff  distance  over  a  small  range  of  Mach  numbers,  while  the  quadratic 
curve  fit  for  the  solution  points  in  Fig.  12  yields  a  similar  dependence  on  the  Mach  number,  AS  [mm]  = 
0.103  M2- 1.151  M  + 4.688. 
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Figure  12.  Standoff  distance  versus  Mach  number  for  various  scenarios  of  air  flow  over  a  6  mm  cylinder. 

Comparing  the  maximum  and  minimum  computed  standoff  distances  to  the  nominal  scenario  shows  the 
uncertainty  in  standoff  distance  due  to  uncertainty  in  freestream  conditions  is  ±3%,  which  is  essentially 
negligible  for  the  energy  deposition  scenario  since  the  maximum  standoff  distance  increase  was  25%. 1  As 
such,  the  nominal  freestream  conditions  listed  in  Table  1  are  employed  for  the  remaining  two-dimensional 
simulations. 


E.  Partial  Slip  Walls 

Since  the  experiment  was  performed  in  a  very  cold,  low  pressure  wind  tunnel,  it  is  possible  that  part  of 
the  cylinder  experiences  velocity  slip  along  its  surface  due  to  the  transition  of  flow  into  to  the  free  molecular 
regime.  To  check  for  partial  velocity  slip,  the  gradient  length  Knudsen  number  developed  by  Boyd  et  al., 
is  computed  from  the  nominal  baseline  steady-state  solution  using  Eq.  (4): 


KnQL  = 


A 

Q 


dQ 

dl 


(4) 
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where  A  is  the  mean  free-path,  Q  is  any  primitive  quantity  (U,  T,  p ,  or  p ),  and  the  derivative  is  taken  in  the 
direction  of  the  maximum  gradient.  Note  when  computing  the  velocity  magnitude  Knudsen  number  gradient, 
the  gradient  is  normalized  by  the  local  velocity  magnitude,  except  when  the  velocity  is  zero.  During  these 
situations,  the  gradient  is  normalized  by  the  local  speed  of  sound.  The  Navier-Stokes  equations  breakdown 
when  the  flow  transitions  from  continuum  to  free  molecular  flow,  which  is  assumed  when  Kiigl  >  0.05. 30 
However,  using  partial  slip  wall  boundary  conditions  where  Kiigl  >  0.05  will  extend  the  range  of  validity  of 
the  Navier-Stokes  solver  being  employed. 

Figure  f3  shows  contours  for  the  gradient  length  Knudsen  number  for  the  nominal  baseline  scenario 
(using  the  grid  shown  in  Fig.  25).  As  seen  in  Fig.  13,  the  flow  experiences  continuum  breakdown  in  the 
vicinity  of  the  shock  and  along  the  top  and  leeward  side  of  the  cylinder.  Continuum  breakdown  within  the 
shock  is  anticipated  due  to  the  strong  flow  discontinuity,  though  breakdown  of  the  governing  equations  in  this 
region  does  not  adversely  affect  the  post-shock  flow  conditions  because  the  Rankine-Hugoniot  condition  is 
preserved  through  the  shock.  However,  continuum  breakdown  observed  on  the  cylinder  surface  is  problematic 
as  it  indicates  that  the  LeMANS  solver  should  not  be  used  in  that  region  without  modifying  the  surface 
boundary  conditions  to  account  for  the  slight  deviation  from  the  continuum  assumption  made  when  deriving 
the  Navier-Stokes  equations  used  in  LeMANS. 


Figure  13.  Contours  of  the  gradient  length  Knudsen  number  for  Mach  5  air  flow  over  a  cylinder. 


In  order  to  quantify  the  effects  of  continuum  breakdown  on  the  shock  structure,  a  partial  slip  wall  model 
is  employed  using  Maxwell’s  slip  condition,  as  described  in  Refs.  30,31.  Maxwell’s  slip  condition  was  derived 
for  a  flat  plate.  It  modifies  the  surface’s  tangential  velocity,  and,  subsequently,  provides  a  temperature  jump. 
Equations  (5)  and  (6)  show  the  modifications  to  the  surface  velocity  and  temperature: 


Uw  =  A 


\  dUx 
d  n 


n 


(5) 


T«,  =  T 


g 


d  T 


(6) 


where  A  is  the  mean  free-path,  T^-T^  is  the  temperature  jump,  a  is  the  accommodation  coefficient,  and  A  is 
a  constant  based  on  the  scenario  being  considered.  For  simplicity  in  this  work,  both  A  and  a  are  set  to  unity, 
as  it  is  unclear  what  they  should  before  hand  for  a  good  nonequilibrium  solution  for  the  flow  considered. 
The  derivative  is  computed  normal  to  the  wall  surface  and  ux  is  the  velocity  in  the  streamwise  direction 
(along  the  wall  surface). 

In  order  to  validate  the  partial  slip  boundary  used,  a  nonequilibrium  solution  is  computed  using  the 
Direct  Simulation  Monte  Carlo  (DSMC),  code:  HAP,  developed  by  Burt  et  al.32  Figure  14  shows  the  surface 
velocity  and  temperature  jump  for  both  the  DSMC  and  LeMANS  with  Maxwell’s  slip  condition.  While 
the  DSMC  solution  exhibits  a  large  amount  of  scatter  (due  to  a  minimal  sample  of  particles  used  in  the 
simulation  and  an  insufficient  number  of  samples  collected  to  further  reduce  the  variance),  the  solution 
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confirms  that  the  flow  experiences  surface  velocity  slip  and  a  small  temperature  jump  as  the  flow  moves  over 
the  top  and  along  the  leeward  side  of  the  cylinder.  In  addition,  the  DSMC  results  provide  direction  for  the 
proper  adjustments  to  A  and  a  in  the  partial  slip  boundary  conditions.  However,  the  values  selected  provide 
reasonable  adjustments  to  the  wall  boundary  conditions  and  are  considered  sufficient  to  quantify  the  effect 
a  partial  slip  wall  has  on  the  bow-shock  location. 


(a)  Tangential  surface  velocity  (b)  Surface  temperature  jump 


Figure  14.  Surface  properties  for  Mach  5  air  flow  over  a  cylinder  from  a  DSMC  simulation  and  LeMANS  using 
a  simple  Maxwell  wall  slip  condition. 


While  the  surface  does  experience  partial  velocity  slip,  which  changes  the  size  and  flow  properties  of 
the  recirculating  wake  region,  the  bow-shock  standoff  is  not  influenced  by  the  change  in  wall  boundary 
conditions,  as  seen  in  Fig.  15.  As  such,  continuum  breakdown  along  the  top  and  leeward  side  of  the  cylinder 
does  not  explain  the  observed  discrepancy  in  bow-shock  standoff  distance  and  is  neglected  in  the  remaining 
simulations. 


Figure  15.  Temperature  contours  for  Mach  5  air  flow  over  a  cylinder  with  and  without  partial  velocity  slip 
and  temperature  jump  wall  boundary  condition. 
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As  seen  in  the  previous  subsections,  uncertainty  in  the  computed  shock  standoff  distance  due  to  freestream 
conditions,  boundary  conditions,  and  thermochemical  effects  is  relatively  small  and  does  not  account  for  the 
discrepancy  in  shock  standoff  seen  in  the  experiment.  While  the  difference  in  shock  standoff  distance  in 
the  baseline  flow  is  important  for  replicating  and  quantifying  the  effects  of  the  discharge  event,  the  next 
section  uses  the  two-dimensional  approach  to  qualitatively  explore  the  nanosecond  DBD  discharge  effect 
and  subsequent  perturbation  on  the  flow.  These  computations  are  performed  to  develop  sufficient  grid  and 
time-step  bounds  for  the  discharge  simulations,  and  to  develop  and  optimize  the  shape  of  the  deposition 
model  before  transitioning  to  three-dimensional  simulations. 

IV.  2D  Energy  Deposition 

The  nanosecond  DBD  event  performed  in  the  experiment  results  in  rapid  thermalization  of  the  electrical 
energy  in  the  surrounding  air,  which  produces  a  compression  wave  that  interacts  with  the  cylinder  bow- 
shock.  Since  the  thermal  efficiency  of  the  DBD  is  not  exactly  known,  one  of  the  critical  parameters  in 
the  reduced-order  model  used  in  the  CFD  modeling  of  event  is  the  rate  of  energy  thermalization,  which  is 
controlled  by  the  parameter  Q  in  Eq.  (1). 

Figure  16  shows  a  grid  used  in  the  energy  deposition  simulation.  The  grid  points  are  spaced  uniformly 
both  in  the  x  and  9  directions  to  reduce  the  spatial  error  due  to  the  propagation  of  the  bow-shock  perturbed 
by  the  compression  wave  formed  during  the  discharge  event.  Three  doubly  refined  grids  have  been  used  to 
study  the  effect  of  grid  resolution  on  the  numerical  solution,  201  x  101  points  (coarse),  401  x  201  points 
(medium)  and  801  x  401  points  (fine).  For  the  coarse  grid,  the  cell  size  is  1  •  10-5  m  x  2.3  •  10-5  m.  The 
computational  domain  length  in  the  streamwise  direction  extends  2  mm  from  the  stagnation  point  of  the 
cylinder  (see  Fig.  16),  which  allows  the  shock  standoff  distance  to  increase  by  up  to  30%  from  the  steady-state 
baseline  solution  standoff  distance  without  being  influenced  by  the  domain  inlet. 

The  spatial  profile  of  energy  deposition  by  the  dis¬ 
charge  used  in  the  calculations  is  set  in  Eq.  (1).  The 
size  of  the  energy  deposition  region  was  determined 
from  the  size  of  the  plasma  visible  emission  in  the  ex¬ 
periment,  and  approximated  as  an  ellipse  with  axes  of 
a  =  0.1  mm  and  b  =  0.3  mm,  with  the  center  of  the 
ellipse  located  at  the  stagnation  point  (for  the  coordi¬ 
nate  system  used  in  the  computations  xc  =  3.0  mm 
and  yc  =  0.0  mm),  and  Gaussian  distribution  of  power 
density.  Thus,  approximately  67%  of  the  input  power 
is  deposited  for  x  <  a  and  y  <  b ,  and  95%  of  power  is 
deposited  for  x  <  2 a  and  y  <  2b.  The  temporal  pro¬ 
file  of  energy  deposition  is  modeled  as  a  step  function. 

Due  to  the  difference  between  baseline  shock  standoff 
distance  predicted  by  the  CFD  model  and  the  experi¬ 
mental  value  (AS  =  1.55  mm  versus  AS  =  1.2  mm),  ob¬ 
taining  agreement  with  all  experimentally  measured  pa¬ 
rameters,  i.e.  compression  wave  speed  of  Uc  =  370  m/s, 
perturbed  bow-shock  propagation  velocity  of  Us  =  92 
m/s,  and  shock  standoff  distance  increase  of  25%,  is 
problematic.  Since  the  compression  wave  is  expand¬ 
ing  radially  as  it  travels,  its  strength  decreases  signif¬ 
icantly  with  distance  traveled  before  interacting  with 
the  bow-shock.  As  a  result,  a  stronger  wave  will  need 
to  be  generated  in  the  computational  exercise  in  order 
to  produce  a  perturbed  bow-shock  with  dynamics  similar  to  those  observed  in  the  experiment. 

To  ensure  spatial  and  temporal  independence  during  the  discharge  event  and  subsequent  flow  interaction, 
simulations  were  conducted  using  an  energy  deposition  pulse  r  =  800  ns  long,  with  total  power  deposited 
Q  =  3  kW,  and  a  fixed  computational  time  step  At  =  0.5  ns.  These  parameters  were  selected  for  the 
grid  independence  study  because  the  deposition  energy  represents  a  high  thermal  efficiency  scenario  and  the 
resultant  flow  perturbation  is  significant. 


X  [m] 


Figure  16.  Grid  used  for  the  energy  deposition  sim¬ 
ulation  (201  x  101). 
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Figure  17(a)  plots  density  gradient  distributions  at  t  =  2  ps  after  the  beginning  of  the  energy  deposition 
pulse,  obtained  using  three  different  grids  and  a  time  step  of  At  =  0.5  ns.  The  density  gradient  was  calculated 
using  a  3-point  central  difference.  The  bow-shock  and  the  propagating  compression  wave  front  are  located 
at  x  =  0.00145  m  and  x  =  0.0021  m,  respectively.  It  can  be  seen  that  the  magnitude  of  the  density  gradient 
increases  approximately  linearly  with  increasing  grid  resolution.  Similarly,  the  magnitude  of  the  peak  density 
gradient  is  also  observed  in  Fig.  17(b),  which  shows  density  gradient  distributions  of  the  perturbed  bow-shock 
at  t  =  6  ps  after  the  beginning  of  the  energy  deposition  pulse.  These  results  indicate  that  the  medium  grid 
has  sufficient  resolution  to  adequately  capture  the  discharge  event  and  resultant  shock-shock  interaction. 


X  [m] 


(b)  t  =  6  [is  (after  shock  interaction) 


Figure  17.  Density  gradient  distributions  for  the  perturbed  bow-shock  obtained  using  three  different  grids  at 
various  times  after  the  energy  deposition  pulse.  Flow  direction  is  left  to  right. 


Figures  18(a)  and  18(b)  show  distributions  of  wall  pressure  and  nondimensional  heat  transfer  coefficients 
on  the  cylinder  model  at  t  =  2  ps  after  the  discharge  pulse,  calculated  for  the  three  different  grids.  It  can 
be  seen  that  with  pulsed  energy  deposition,  peak  wall  pressure  is  reached  off  stagnation  line,  approximately 
at  9  =  20°  (see  Fig.  18(a)).  Pulse  energy  deposition  also  results  in  significant  heat  transfer  increase  near 
the  stagnation  line  (see  Fig.  18(b)).  These  results  indicate  that  all  three  grids  provide  sufficient  resolution 
to  obtain  accurate  estimates  of  the  cylinder  surface  properties. 


(b)  Nondimensional  heat  transfer  coefficient 


Figure  18.  Surface  distributions  at  t  =  2  \xs  after  energy  deposition  pulse  for  various  grids. 
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Figure  19  plots  compression  wave  speed  and  perturbed  bow-shock  velocity  from  simulations  for  three 
fixed  time  steps:  At  =  5  ns,  1  ns,  and  0.5  ns  on  the  three  different  grids.  To  calculate  the  compression  wave 
speed,  the  location  of  the  density  peak  was  sampled  every  400  ns  at  t  =  1  —  3  ps  after  the  beginning  of  the 
energy  deposition  pulse.  After  the  compression  wave  reached  the  standing  bow-shock  (approximately  4  p.s 
after  the  beginning  of  the  energy  deposition  pulse),  the  location  of  the  density  peak  was  sampled  every  400 
ns  at  t  =  4  —  6  (is  after  the  beginning  of  the  energy  deposition  pulse  to  determine  the  perturbed  bow-shock 
speed.  From  Fig.  19,  it  can  be  seen  that  the  calculated  compression  wave  speed  and  perturbed  bow-shock 
speed  depend  on  grid  resolution  when  the  computation  time  step  is  At  =  5  ns,  but  fully  converges  when 
the  time  step  is  reduced  to  At  =  0.5  ns.  Based  on  these  results,  a  fixed  time  step  of  At  =  0.5  ns  and  the 
medium  grid  (401  x  201),  provide  sufficient  temporal  and  spatial  independence  and  is  used  for  all  subsequent 
two-dimensional  simulations. 

Figure  20  shows  dependence  of  compression  wave  propagation  velocity  on  deposited  pulse  energy  for 
input  powers  of  Q  =  1,  2,  and  3  kW.  Since  the  energy  is  deposited  uniformly  in  time,  the  duration  of  energy 
deposition  pulse  is  varied  as  follows:  400ns  <  r  <  3200  ns  for  Q  =  1  kW,  200ns  <  r  <  1600  ns  for  Q  =  2 
kW,  and  200ns  <  r  <  1200  ns  for  Q  =  3  kW,  respectively.  Fig.  20  also  denotes  the  speed  of  sound  evaluated 
using  the  post-shock  temperature  averaged  along  the  stagnation  line  (T=  307  K,  Cs  =  351  m/s). 


Figure  19.  Compression  wave  propagation  veloc-  Figure  20.  Compression  wave  speed  dependence 

ity  (closed  symbols)  and  perturbed  bow-shock  ve-  on  input  power  and  energy  added  (i.e.  duration 

locity  (open  symbols)  vs.  time  step  for  the  three  of  energy  deposition  pulse), 

grids  used. 

As  seen  in  Fig.  20,  input  power  of  Q  =1  kW  and  2  kW  produce  a  compression  wave  speed  that  increases 
with  total  pulse  energy  until  an  energy  of  1.6  mJ,  whereas  the  input  power  of  Q  =  3  kW  achieves  a 
compression  wave  speed  that  increases  up  to  pulse  energy  2.4  mJ.  This  is  due  to  the  fact  that  the  smaller 
values  of  Q  fail  to  produce  a  supersonic  compression  wave,  which  limits  the  amount  of  energy  that  can 
be  used  to  increase  the  wave  speed.  Since  the  compression  wave  is  supersonic  in  the  experiment  (see  Fig. 
21),  a  value  of  Q  >  3  kW  seems  most  appropriate  to  investigate  the  discharge  effect  and  subsequent  flow 
perturbation.  In  addition,  selection  of  an  input  power  of  Q  =  3  kW  for  r  =  800  ns,  corresponds  to  a 
pulse  energy  of  2.4  mJ,  which  is  consistent  with  about  30%  the  experiment’s  coupled  pulse  energy  ‘quickly’ 
thermalizing.  The  thermal  efficiency  assumed  here  is  consistent  with  a  study  by  Popov.22 

Figure  22(a)  shows  density  distributions  along  the  stagnation  streamline,  calculated  for  pulse  duration 
of  r  =  1200  ns  (Q  =  3  kW),  at  t  =  200  ns,  400  ns,  800  ns,  and  1200  ns  after  the  beginning  of  the  energy 
deposition.  As  seen  in  the  figure,  the  region  of  fluid  displacement  (i.e.  density  reduction  near  the  stagnation 
point)  expands  up  to  t  =  800  ns,  which  results  in  density  increase  to  the  left  of  the  energy  deposition 
region.  At  t  >  800  ns,  the  density  in  the  energy  deposition  region  no  longer  decreases,  thus  limiting  mass 
displacement  due  to  energy  deposition  and,  consequently,  compression  wave  speed.  After  this  limit  is  reached, 
additional  energy  input  no  longer  contributes  to  increasing  the  compression  wave  speed.  As  a  result,  the 
compression  wave  speed  remains  constant  for  pulse  deposition  durations  greater  than  800  ns  (for  Q  =  3  kW) . 
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Figure  21.  Comparison  of  predicted  peak  density  gradient  location  with  the  experiment  (phase- locked 
schlieren).  (401  x  201,  At  =  0.5  ns,  Q  =  3  kW,  r  =  800  ns) 


X  [m] 


(a)  Density 


Figure  22.  Stagnation  streamline  distributions  at  different  delays  after  the  start  of  the  energy  deposition  pulse 
(Q  =  3  kW,  t  =  800  ns),  for  Mach  5  air  over  a  cylinder. 


Using  a  total  input  power  of  Q  =  3  kW  and  energy  deposition  pulse  duration  of  r  =  800  ns  yields  the 
maximum  compression  wave  speed  for  its  total  input  power  while  providing  a  wave  speed  that  is  slightly 
greater  than  the  experiment’s  measured  wave  speed.  As  such,  this  scenario  provides  the  best  agreement 
for  both  compression  wave  velocity  and  perturbed  bow-shock  propagation  velocity  along  the  stagnation 
line,  as  shown  in  Fig.  21.  The  total  energy  deposited  for  this  scenario  amounts  to  2.4  mJ/pulse  (for  the 
full  geometry),  which  means  the  nanosecond  DBD  pulse  would  have  to  have  a  thermal  efficiency  of  over 
30%.  However  this  conclusion  can  be  confusing  since  the  baseline  bow-shock  location  is  1.55  mm  from  the 
cylinder  in  the  computation  and  only  1.2  mm  for  the  experiment.  Since  the  shock  is  located  farther  from 
the  cylinder,  the  energy  deposited  in  the  simulation  must  be  exaggerated  in  order  to  create  a  compression 
wave  with  sufficient  speed  and  strength  such  that  its  interaction  with  the  standing  bow-shock  produces  a 
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resultant  combined  wave  that  has  the  same  velocity  as  the  experiment. 

Figure  22(b)  plots  temperature  distributions  at  the  same  delays  as  Fig.  22(a).  In  these  calculations, 
temperature  at  the  stagnation  point  is  fixed  at  T  =  300  K  (isothermal  wall).  However,  the  peak  temperature 
increases  during  the  energy  deposition  process  and  reaches  maximum  at  t  =  800  ns,  T  =  2160  K.  This 
peak  temperature  value  is  much  higher  than  the  temperature  measured  in  the  experiment  using  N2  second 
positive  bands  emission  spectra,  T  =  340  ±  30  K.  However,  it  is  well  known  that  N2  second  positive  system 
emission  decays  over  a  few  tens  of  nanoseconds  after  the  discharge  pulse.  Note  that  temperature  in  the  energy 
deposition  region  predicted  at  t  =  40  ns  more  reasonably  matches  the  experiment,  with  a  peak  temperature 
of  T  =  410  K.  In  addition,  the  average  temperature  along  the  stagnation  line  (within  the  shock  layer)  during 
the  first  40  ns  of  the  energy  deposition  pulse  computation  is  T  =  355  K. 

Confirmation  of  the  post-shock  temperature  50  nanoseconds  after  the  discharge  event  suggests  that  the 
significant  temperature  rise  seen  later  in  the  computation  could  actually  be  occurring  in  the  experiment  or 
that  the  strong  temperature  rise  is  merely  an  artifact  of  the  phenomenological  deposition  modeling  processing 
including  the  shape,  distribution,  and  duration  of  the  deposition.  Verifying  the  this  requires  taking  additional 
temperature  measurements  during  the  experiment  at  much  longer  time  delays  after  the  discharge  pulse  (t  =  1 
ps),  a  challenge  still  being  pursued,  though  the  recent  picosecond  CARS  results33  may  produce  a  method 
for  obtaining  the  necessary  measurements. 

Figure  23(a)  plots  wall  pressure  coefficient  distributions  for  several  time  delays  after  the  energy  deposition 
pulse.  It  can  be  seen  that  the  initial  pressure  rise  due  to  energy  deposition  is  followed  by  a  rarefaction  wave. 
At  t  =  1.2  ps  and  2.0  ps,  the  rarefaction  wave  reduces  the  wall  pressure  coefficient  below  baseline  value 
for  the  first  15°  along  the  cylinder  surface.  The  wall  pressure  perturbation  is  subsequently  reduced  and 
approaches  baseline  profile  approximately  at  t  =  4  ps,  which  has  been  excluded  from  figure  as  it  overlaps 
the  baseline  solution. 


Figure  23.  Surface  pressure  coefficient  for  different  time  delays  after  beginning  of  energy  deposition  pulse 


Figure  23(b)  shows  wall  pressure  coefficient  distributions  at  t  =  8  ps  and  9  ps,  displaying  a  second  wall 
pressure  ‘dip’  near  the  stagnation  line.  The  second  pressure  reduction  is  due  to  a  rarefaction  wave,  which 
reflects  off  the  perturbed  bow-shock  and  returns  back  to  the  cylinder.  The  movement  of  the  rarefaction  wave 
is  also  visible  in  Fig.  24,  which  plots  pressure  distributions  along  the  stagnation  streamline  for  different  time 
delays  after  the  pulse. 

To  determine  the  drag  on  the  cylinder,  the  computational  domain  is  extended  to  include  the  recirculating 
wake  region,  as  shown  in  Fig.  25.  The  grid  used  for  these  calculations  is  401  x  401  points,  with  the  grid 
resolution  for  0°  <  6  <  90°  being  the  same  as  the  medium  resolution  grid  discussed  previously. 

Figure  26(a)  shows  time  evolution  of  total  drag  (i.e. ,  the  sum  of  pressure  drag  and  shear  drag  over  the 
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Figure  24.  Pressure  distributions  along  the  stagnation  streamline  at  t  =  5  |as,  7  |j.s,  and  8  jj.s  delays  after  the 
energy  deposition  pulse. 


cylinder  surface).  The  total  drag  for  the  baseline  steady-state  scenario  is  D  =  2.4  N.  Note  that  shear  drag 
accounts  for  1%  of  the  drag  due  to  the  low  temperatures  observed  in  the  experiment.  The  drag  reaches  a 
maximum  D  =  3.04  N  (26%  increase  from  the  baseline),  at  t  =  400  ns,  and  a  minimum  D  =  2.27  N  (5% 
decrease  from  the  baseline)  at  t  =  8.5  ps.  Following  the  interaction  with  the  rarefaction  wave  on  the  cylinder, 
the  drag  recovers  the  baseline  value  around  t  =  15  ps  after  the  discharge  event. 
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Figure  25.  Computational  domain  used  for  drag  calculations  (401  X  401). 


To  evaluate  the  effect  of  energy  deposition  on  drag  reduction,  the  time-averaged  change  in  total  drag  is 
shown  in  Fig.  26(b).  As  seen  in  the  figure,  the  drag  increases  up  to  20%  at  t  =  1  ps,  and  then  monotonically 
decreases  to  -1%  at  t  =  11  ps,  and  remains  unchanged  thereafter.  Figure  26(b)  show  that  the  nanosecond 
DBD  discharge  provides  a  slight  reduction  in  total  drag  over  the  lifetime  of  the  discharge  event.  However, 
the  drag  reduction  comes  with  a  substantial  increase  to  the  heat  load  on  the  cylinder. 

To  test  if  there  is  a  continued  reduction  in  drag  by  repetitive  discharge  events,  three  additional  scenarios 
are  considered.  The  first  repeats  the  energy  deposition  at  t  =  15  ps  (i.e,  the  discharge  event  begins  after 
the  previous  discharge  event  cycle).  This  would  be  consistent  with  an  operational  frequency  of  67  kHz.  The 
results  essentially  replicate  those  previously  presented  for  a  single  pulse  effect,  which  suggests  the  1%  drag 
reduction  on  the  cylinder  could  be  sustained  if  the  DBD  was  pulsed  at  this  frequency.  In  the  second  scenario, 
the  second  discharge  event  begins  when  the  total  drag  is  at  its  lowest  (t  =  8.5  ps  after  the  discharge  event), 
which  corresponds  to  an  operational  frequency  of  110  kHz.  The  third  scenario  starts  the  second  discharge 
event  at  t  =  11  ps,  which  corresponds  to  the  lowest  time-averaged  drag,  or  an  operating  frequency  of  90 
kHz. 

Figure  27  shows  the  total  drag  time  histories  for  repetitive  discharges.  As  seen  in  the  figure,  neither 
of  the  higher  frequency  repetitive  energy  deposition  scenarios  resulted  in  an  improvement  to  the  total  drag 
reduction.  This  occurs  because  the  temperature  in  the  energy  deposition  region  is  slightly  hotter  during 
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(b)  Percent  change  of  the  ‘running  average’  drag 
obtained  by  time  integration  of  Fig.  26(a)  nor¬ 
malized  by  the  baseline  flow. 


Figure  26.  Time  evolution  of  drag  on  the  cylinder,  calculated  using  the  grid  show  in  Fig.  25. 


the  second  discharge  (for  the  two  higher  frequency  scenarios),  since  the  flow  has  less  time  to  cool  down 
between  discharge  events.  As  a  result,  the  energy  deposited  during  the  second  pulse  has  a  diminished  effect 
in  generating  the  subsequent  compression  wave,  so  its  interaction  with  bow-shock  is  also  reduced.  These 
results  suggest  that  it  is  best  to  operate  the  nanosecond  DBD  actuator  at  67  kHz  when  it  is  being  used  for 
continuous  operation  since  the  total  power  requirements  are  the  lowest,  and  it  is  still  able  to  maintain  the 
1%  reduction  in  running  average  total  drag  over  the  discharge  cycle. 


Figure  27.  Time  evolution  of  the  total  drag  on  the  cylinder,  with  a  second  energy  deposition  pulse  generated 
at  t  =  9ps  (110  kHz)  and  t  =  11  ps  (90  kHz). 


As  a  final  exploration  using  two-dimensional  simulations,  a  simulation  is  performed  with  an  input  power 
Q  =  35  kW  over  r  =  200  ns  (7  mJ/pulse  for  a  whole  cylinder,  though  the  computation  was  performed 
using  the  1/a  cylinder  domain).  This  simulation  was  conducted  in  order  to  obtain  a  compression  wave  with 
sufficient  speed  such  that  the  shock-shock  interaction  occurred  in  the  same  time  as  the  experiment.  In 
Fig.  28,  simulation  results  are  shown  at  the  same  time  moments  as  the  experimental  phase-locked  schlieren 
images.  The  figure  shows  many  similarities  between  the  bow-shock  perturbation  from  the  experiment  and 
the  computation. 
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(a)  Phase-locked  schlieren  images 


(b)  CFD  generated  density  gradient  distributions 


Figure  28.  Temporal  evolutions  of  the  bow-shock  (a)  schlieren  photography  (experiment)  and  (b)  density 
gradient  distribution  in  the  streamwise  direction  (CFD). 


As  seen  in  the  Fig.  28,  the  generated  compression  wave  propagates  upstream  and  reaches  the  bow-shock 
at  t  =  3  ps.  The  interaction  results  in  a  25%  increase  in  the  shock  location  from  the  baseline  standoff 
distance,  which  is  consistent  with  the  results  observed  in  the  experiment’s  schlieren  images.  Note  in  this 
simulation,  higher  power  input  was  required  for  the  compression  wave  to  reach  the  bow-shock  at  t  =  3  ps, 
while  compensating  for  a  30%  difference  in  the  bow-shock  standoff  distance.  It  would  require  a  significantly 
smaller  amount  of  deposition  energy  for  a  compression  wave  to  travel  1.2  mm  in  the  same  amount  of  time. 

V.  3D  Simulations 

As  previously  discussed,  the  two-dimensional  nominal  baseline 
simulation  predicts  a  bow-shock  standoff  distance  that  is  signifi¬ 
cantly  larger  than  that  observed  in  the  experiment.  Several  two- 
dimensional  simulations  were  performed  to  quantify  the  uncertainty 
of  the  freestream  and  boundary  conditions,  along  with  thermo¬ 
chemical  nonequilibrium  flow  effects.  However,  the  results  did  not 
explain  the  observed  discrepancy  in  shock  location.  As  such,  a  sim¬ 
ulation  of  the  entire  tunnel  (including  the  nozzle  and  the  region 
downstream  of  the  test  section),  is  performed  in  order  to  replicate 
the  experiment’s  shock  standoff  distance  and  identify  the  reason  for 
the  discrepancy  in  bow-shock  location.  Details  about  the  tunnel  were 
previously  mentioned  in  Section  II.  Note  Fig.  29  provides  an  outline 
of  the  computational  domain  used  in  the  simulation.  Only  1/i  of  the 
tunnel  is  simulated  because  the  tunnel  is  assumed  to  be  symmetric 
in  both  the  spanwise  and  transverse  directions.  The  nozzle  throat  conditions  used  in  the  simulation  are  listed 
in  Table  5. 


Table  5.  Input  conditions  for  the  Mach 
5  wind  tunnel  at  the  nozzle  throat. 


Parameter 

Value 

u* ,  [m/s] 

318 

T*,  [K] 

250 

P *,  [kg/m3] 

0.3609 

p*,  [Pa] 

26,020 

Mach* 

1.0 
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Figure  29.  The  computational  domain  for  the  entire  tunnel.  The  gray  surface  depicts  the  tunnel  while  the 
pink  surface  illustrates  the  location  of  the  cylinder  test  article.  Only  1/4  of  the  geometry  is  simulated. 


Using  the  lessons  learned  from  the  two-dimensional  simulations,  a  grid  was  developed  for  the  tunnel 
using  4  computational  blocks.  Each  block  is  composed  of  structured  cells.  The  grid  spacing  is  such  that  grid 
clustering  occurs  near  all  surfaces.  The  computational  block  surrounding  the  cylinder  geometry  is  identical 
to  the  grid-independent  ‘medium’  grid  used  for  the  two-dimensional  analysis.  In  total,  the  three-dimensional 
computational  domain  contains  15  M  cells  and  was  run  using  512  processors.  For  the  baseline  simulation 
(i.e.,  without  the  discharge  event),  implicit  time  integration  was  employed  with  a  time  step-size  varying  from 
At  =  0.1  ns  to  At  =  2  qs. 

No  slip,  isothermal  walls  (T^  =  300  K),  are  assumed  for  all  surfaces  and  a  non-reflective  first-order 
extrapolation  is  used  at  the  domain  exit  plane.  The  simulation  is  started  from  quiescent  air,  except  for  the 
input  conditions  listed  in  Table  5.  As  a  result,  the  simulation  required  about  20  milliseconds  (~  30, 000 
iterations),  for  the  wall  boundary  layer  and  cylinder  bow-shock  structure  to  develop  and  achieve  a  quasi 
steady-state  in  the  inviscid  region.  Figure  30  shows  a  Mach  1.25  iso-surface  of  the  flow,  which  is  colored 
by  temperature.  This  iso-surface  level  was  selected  as  it  illustrates  the  substantial  growth  of  the  sidewall 
boundary  layer  as  the  flow  expands  in  the  tunnel.  It  is  also  easy  to  see  that  the  top  wall  boundary  does  not 
grow  as  dramatically  due  to  the  1.5°  divergence  of  the  top  wall.  In  addition,  the  figure  shows  the  bow-shock 
standoff  distance  is  1.23  mm  from  the  cylinder  (as  seen  in  the  zoomed  box),  which  matches  the  experiment’s 
measurement  of  the  shock  location. 


Figure  30.  Mach  1.25  iso-contour  colored  by  temperature  at  t  =  19  ms.  The  pink  surface  indicates  the  cylinder 
for  the  1/4  geometry  simulation  of  air  expanding  in  a  Mach  5  wind  tunnel. 
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Figure  31  shows  the  phase-locked  schlieren  images  from  the  experiment,  which  have  been  overlaid  with 
the  computational  solutions.  The  schlieren  images  in  the  experiment  were  taken  with  a  knife-edge  set  in  the 
streamwise  direction.  To  replicate  this,  the  computational  schlieren  result  is  computed  in  the  same  direction 
(i.e.,  spanwise  average  of  the  derivative  of  density  in  the  x-direction) .  Accounting  for  the  thickness  and 
location  of  the  Kapton®  tape  placed  on  the  cylinder,  it  can  be  seen  that  the  computational  shock  standoff 
distance  matches  the  experimental  quite  well. 


(a)  Side  view 


(b)  Top-down  view 


Figure  31.  Schlieren  images  from  the  experiment,  along  with  the  computational  density  gradient  in  the  x- 
direction  at  t  =  40  ms. 


Aside  from  schlieren  visualization,  the  only  other  validation  available  from  the  experiment  is  the  surface 
pressure  tap  located  on  sidewall  of  the  tunnel,  4  cm  upstream  of  the  cylinder.  Figure  32  plots  the  pressure 
on  the  tunnel  sidewall.  The  computed  pressure  at  the  location  of  the  pressure  tap  is  pw  =  1.1  Torr,  which  is 
slightly  lower  than  the  experimental  value  of  pw  =  1.2  Torr.  However,  given  the  precision  of  the  experimental 
measurement  (±0.05  Torr),  and  the  uncertainties  associated  with  the  tunnel  conditions,  the  computational 
results  are  consider  to  be  in  very  good  agreement  with  the  experimental  measurement. 


0.04 


^0.02 

> 


X  [m] 


Figure  32.  Pressure  contours  along  the  sidewall  the  Mach  5  wind  tunnel  with  a  6  mm  cylinder  in  the  test 
section.  Flow  direction  is  left  to  right.  ( t  =  20  ms) 
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Iii  addition  to  the  location  of  the  bow-shock,  Fig.  30  also  illustrates  the  complexity  of  the  flow-field, 
including  a  reflected  Mach  wave  downstream  of  the  cylinder  wake  along  the  flow  spanwise  centerline.  In 
particular,  the  interaction  between  the  cylinder  and  sidewall  results  in  a  complex  wake  region  that  oscillates 
at  a  low  frequency  (30  Hz).  The  low  frequency  oscillation  of  the  boundary  layer  is  due  to  interaction  between 
the  sidewall  and  the  cylinder  body,  and  is  a  characteristic  of  wall/blunt-body  experiments,  where  a  lambda 
shock  structure  forms  at  the  edge  of  the  inviscid  region  and  the  boundary  layer.  The  oscillation  in  the  wake 
region  of  the  cylinder /side- wall  interaction  can  be  seen  in  Fig.  33,  which  shows  the  Mach  1.25  iso-surface,  but 
viewed  looking  upstream.  Oscillations  at  the  shock/boundary  layer  junction  were  also  observed  in  schlieren 
images  from  the  experiment. 

As  seen  in  Fig.  34,  the  pressure  is  essentially  uniform  in  the  boundary  layer,  but  rapidly  changes  in  the 
inviscid  region  (decreasing  when  the  Mach  number  is  increasing  and  increasing  where  the  Mach  number 
is  decreasing).  As  such,  the  assumptions  made  to  originally  estimate  the  tunnel’s  flow  conditions  were 
appropriate  (i.e.,  constant  static  pressure  through  the  boundary  layer),  but  the  approach  did  not  account  for 
the  rapid  increase  in  boundary  layer  thickness  so  far  upstream  of  the  test  article.  As  a  result,  the  inviscid  flow 
just  upstream  of  the  cylinder’s  bow-shock  is  Mach  4  (rather  than  Mach  5).  According  to  Billig’s  empirical 
formulation,  the  decrease  in  freestream  Mach  number  should  further  increase  the  shock  standoff  distance, 
but  due  to  the  complex  flow  structures  from  the  cylinder/side- wall  junction,  the  bow-shock  remains  close  to 
the  geometry. 


(c)  t  =  30  ms 


(d)  t  =  40  ms 


Figure  33.  Mach  1.25  iso-surface  colored  by  temperature  for  Mach  5  air  in  a  wind  tunnel.  The  cylinder /side- 
wall  junction  results  in  a  complex  wake  flow  that  ‘breathes’  over  a  30  ms  cycle  (~  30  Hz). 


As  the  flow  accelerates  through  the  nozzle,  the  pressure  in  the  inviscid  core  drops,  but  due  to  the  thick 
boundary  layer  far  upstream  of  the  test  section,  the  inviscid  core  flow  eventually  contracts  and  slows  slightly. 
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The  decrease  in  Mach  number  coincides  with  an  increase  in  static  pressure  within  the  inviscid  core.  As  a 
result,  there  is  a  higher  post  bow-shock  pressure  and  a  larger  spanwise  pressure  gradient  in  that  region. 
Figure  34  shows  Mach  number  and  pressure  contours  along  the  transverse  symmetry  plane  upstream  of  the 
test  section. 
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Figure  34.  Top-down  view  of  Mach  and  pressure  contours  along  the  transverse  symmetry  plane  (y  =  0),  for 
air  flow  in  a  Mach  5  wind  tunnel. 


Figure  35  shows  a  slice  of  pressure  contours  along  the  transverse  symmetry  plane  near  the  cylinder.  As 
seen  in  the  figure,  there  is  a  very  large  pressure  drop  in  the  spanwise  direction  just  after  the  bow-shock.  This 
pressure  drop  is  due  to  the  low  pressure  present  in  the  boundary  layer,  and  is  further  strengthened  due  to 
the  strong  interaction  between  the  cylinder  with  the  sidewall. 


Figure  35.  Top-down  view  of  pressure  contours  along  the  transverse  symmetry  plane  (y  =  0),  for  air  flow  in 
a  Mach  5  wind  tunnel.  The  figure  also  includes  streamlines  of  velocity.  Note  contours  p  <  1500  Pa  have  been 
excluded. 


The  large  pressure  gradient  causes  the  inviscid  core  flow  to  turn  into  the  spanwise  direction  after  it  goes 
through  the  bow-shock.  While  the  post-shock  flow  escaping  into  the  boundary  layer  is  subsonic,  it  still  has 
a  large  velocity  and,  as  such,  a  large  portion  of  the  mass  flow  is  no  longer  two-dimensional.  As  a  result  of 
the  three-dimensionality  of  the  flow,  the  bow-shock  location  is  much  closer  to  the  cylinder  than  would  be 
expected  in  a  two-dimensional  flow. 
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A.  Energy  Deposition 

Once  the  three-dimensional  baseline  simulation  results  were  completed,  as  simulation  of  the  discharge  event 
was  computed  using  the  same  phenomenological  discharge  model  as  the  two-dimensional  results,  except  the 
model  was  extended  for  three-dimensions: 


x  =  (x-  xc)  (7) 

y  =  (y-  2/c) 
z  =  (z-  zc) 

where  the  center  of  the  ellipsoid  is  set  to  the  stagnation  point  of  the  cylinder  along  both  symmetry  planes  (so 
half  the  ellipsoidal  volume  is  inside  the  cylinder,  and  is  excluded  from  the  computational  domain).  For  the 
wind  tunnel  geometry  shown  previously,  the  values  are  x  =  14.2  cm,  y  =  0,  z  =  2.0  cm.  The  equatorial  radii 
of  the  ellipsoidal  deposition  are  consistent  with  the  two-dimensional  simulations,  and  the  polar  radius,  c,  is 
half  the  width  of  the  exposed  electrode  (a  =  0.1  mm,  b  =  0.3  mm,  c  =  5  mm).  Note  f f f_QO  S  dx  dy  dz  =  Q 
The  simulation  was  carried  out  assuming  a  total  power  deposition  Q  =  5  kW  for  100  ns  (for  1/i  the 
geometry),  at  a  time  step  of  At  =  0.5  ns  to  ensure  temporal  independence  based  on  the  two-dimensional 
simulations.  The  power  deposited  amounts  to  a  total  energy  deposition  of  2  mJ/pulse,  (i.e.,  a  thermal 
efficiency  of  over  25%).  While  the  computational  representation  of  the  discharge  event  results  in  a  relatively 
high  thermal  efficiency,  the  value  of  Q  was  selected  to  ensure  a  supersonic  compression  wave  was  generated. 

The  simulation  is  run  for  15  ps  using  a  At  =  0.5  ns  time-step  (i.e.,  30,000  iterations)  to  capture  the 
evolution  of  the  compression  wave/bow-shock  interaction.  Figure  36  shows  the  Mach  1.25  iso-surfaces  colored 
by  temperature  at  t  =  4.1  ps  after  the  deposition  event. 


Figure  36.  Mach  1.25  iso-contour  colored  by  temperature  for  air  in  a  Mach  5  wind  tunnel  t  =  4.1  ps  after  a 
discharge  event.  The  compression  wave  pushes  the  bow-shock  outward,  as  seen  in  the  red  region. 


Consistent  with  the  two-dimensional  results,  the  power  deposited  by  the  phenomenological  model  was 
sufficiently  high  to  generate  a  supersonic  compression,  which  traveled  upstream  and  interacts  with  the 
standing  bow-shock.  However,  either  the  energy  was  deposited  over  too  small  an  area,  and/or  too  much 
energy  was  deposited,  because  the  compression  wave  speed  was  slightly  faster  than  the  experimental  value 
and,  thus,  started  interacting  with  the  bow-shock  before  the  time  observed  experimentally  with  the  phase- 
locked  schlieren  images.  Further  evolution  of  the  flow  shows  the  resultant  compression  wave/bow-sliock 
structure  extending  out  more  than  2.2  mm  from  the  cylinder  (along  the  centerline),  which  is  a  85%  increase 
in  the  shock  standoff  distance.  This  increase  is  much  higher  than  that  observed  in  the  experiment  (25%).  In 
addition,  the  shape  of  the  perturbed  shock  does  not  completely  match  the  shape  in  schlieren  images.  These 
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results  suggest  that  both  the  shape  and  amount  of  energy  deposited  in  the  phenomenological  deposition 
model  need  to  be  adjusted  in  order  to  fully  replicate  the  experiment  using  three-dimensional  simulations. 

The  extra  computations  needed  to  fine  tune  the  three-dimensional  phenomenological  energy  deposition 
model  are  beyond  the  current  resources  dedicated  for  this  work.  Future  studies  will  further  refine  the  reduced- 
order  model,  though  given  the  high  computational  cost  associated  with  the  three-dimensional  simulations, 
it  may  be  more  feasible  to  experimentally  explore  the  nanosecond  DBD  on  a  5  mm  diameter  cylinder  since 
the  two-dimensional  computations  and  two-dimensional  empirical  correlation  are  in  much  better  agreement 
with  the  experiment. 

While  the  three-dimensional  computation  results  do  not  replicate  the  exact  behavior  observed  in  the 
experiment,  the  solutions  do  show  that  the  compression  wave’s  interaction  with  the  bow-slrock  is  directly 
responsible  for  the  bow-shock  movement  and  that  the  movement  is  not  the  result  of  a  portion  of  the  compres¬ 
sion  wave  interacting  the  with  flow  in  the  cylinder/side- wall  junction,  since  the  flow  remains  quasi-steady  in 
that  region  during  the  entire  15  ps  discharge  cycle. 

VI.  Conclusions 

A  computational  study  of  Mach  5  air  flow  around  a  6  mm  cylinder  with  a  nanosecond  DBD  discharge 
was  replicated  using  high-fidelity  numerical  simulations  in  order  to  better  understand  the  resulting  flow 
perturbation  and  the  impact  the  discharge  event  has  on  the  flow  and  the  cylinder  surface.  The  shock 
standoff  distance  for  the  6  mm  cylinder  scenario  was  found  to  be  20%  smaller  in  the  experiments  than  the 
predictions  of  the  empirical  correlation  and  computations.  However,  the  computed  shock  standoff  distance 
for  a  5  mm  cylinder  agreed  well  with  the  empirical  correlation  and  experiment.  In  order  to  determine  the 
cause  of  the  shock  location  discrepancy,  the  tunnel’s  uncertainty  in  freestream  conditions  were  estimated, 
and  computations  performed.  The  freestream  differences  resulted  in  a  variation  of  the  bow-shock  standoff 
distance  of  about  ±3%.  In  addition,  the  effects  due  to  wall  boundary  conditions,  including  partial  slip  walls, 
as  well  as  thermo-chemical  nonequilibrium  were  found  to  provide  an  insignificant  change  in  the  bow-shock 
location.  Therefore,  the  observed  discrepancy  is  attributed  to  the  interaction  at  the  junction  of  the  cylinder 
and  sidewall  which  results  in  a  complex,  quasi-steady  boundary  layer.  This  complex  flow  acts  to  lower  the 
post-shock  pressure  thereby  drawing  the  bow-shock  closer  to  the  cylinder.  This  was  verified  with  a  three- 
dimensional  simulation  of  the  entire  wind  tunnel,  which  was  able  to  computationally  replicate  the  bow-slrock 
structure  seen  in  the  schlieren  photography,  predict  the  width  of  the  tunnel’s  inviscid  core,  and  match  static 
pressure  with  the  experiment’s  sidewall  pressure  tap. 

The  effect  of  the  nanosecond  pulse  surface  dielectric  barrier  discharge  on  a  Mach  5  flow  over  a  cylinder 
model  was  modeled  using  two-dimensional  compressible  Navier-Stokes  flow  code,  LeMANS,  which  incorpo¬ 
rated  a  phenomenological  thermal  energy  deposition  model.  In  spite  of  a  difference  between  the  baseline 
bow-shock  standoff  distance  measured  in  the  experiments  and  that  predicted  by  the  CFD  model,  other  ex¬ 
perimentally  measured  parameters  such  as  compression  wave  speed  and  perturbed  bow-shock  speed,  as  well 
as  the  time  evolution  of  the  perturbed  bow-shock  shape  were  replicated  using  the  model.  Distributions  of 
wall  pressure  and  heat  transfer  coefficient,  as  well  as  total  drag  on  the  model  have  also  been  calculated  and 
show  that  the  cylinder  experiences  a  small  reduction  in  total  drag  (1%),  over  a  discharge  event  cycle  of  15 
ps,  though  the  cylinder  is  subjected  to  a  large  increase  in  peak  heating  during  the  first  microsecond  of  the 
discharge  cycle. 

A  three-dimensional  simulation  was  also  computed  to  replicate  the  discharge  event  from  the  baseline  flow 
conditions.  While  these  results  do  not  agree  as  well  with  the  experiment’s  phase-locked  schlieren  images 
as  the  two-dimensional  results,  the  resultant  flow  field  indicates  that  the  cylinder/sidewall  junction  is  not 
immediately  influenced  by  the  discharge  event  and  therefore  is  not  responsible  for  the  increase  in  bow-shock 
standoff  distance  due  to  the  nanosecond  DBD  pulse.  In  addition,  the  behavior  of  the  shock  perturbation 
is  independent  of  the  form  of  energy  input  (due  to  the  short  pulse  duration),  so  the  resultant  flow  is  being 
accurately  simulated. 

While  the  nanosecond  DBD  actuator  simulated  in  this  work  is  not  very  effective  at  improving  the  cylinder 
surface  conditions  (i.e.,  it  creates  a  large  increase  in  heat  load  on  the  cylinder  surface  with  only  a  minor 
decrease  to  total  drag),  the  device  is  very  effective  at  moving  a  strong  standing  shock.  This  technology  could 
prove  very  useful  in  supersonic  inlets  and  isolators  found  in  RAM  jet  and  SCRAM  jet  engines,  where  engine 
‘un-start’  continues  to  be  a  technical  challenge.  One  example  of  engine  ‘un-start’  is  the  upstream  movement 
of  the  shock-train  within  the  inlet  which  leads  to  unfavorable  conditions  in  the  combustor.  If  ‘un-start’  of 
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the  shock-train  could  be  sensed  in  the  inlet,  a  nanosecond  DBD  device  could  be  rapidly  fired  to  move  the 
shock-train  back  to  its  nominal  location. 

From  this  work,  it  appears  these  devices  would 
be  best  situated  in  subsonic  regions  of  the  flow  (i.e., 
in  the  sub-sonic  region  of  a  separation  bubble  or 
in  the  subsonic  region  downstream  of  a  bow-shock), 
where  a  supersonic  compression  wave  could  more 
easily  be  generated.  Figure  37  illustrates  the  con¬ 
cept  of  placing  a  nanosecond  DBD  actuator  at  the 
leading  edge  of  a  supersonic  inlet  flow-path.  As  seen 
in  the  figure,  the  DBD  induced  bow-shock  pertur¬ 
bation  travels  along  the  standing  shock,  effectively 
moving  the  subsequent  shock-train. 

This  concept  already  has  received  some  atten¬ 
tion.  In  2008,  Gnemmi  et  al.  successfully  applied  the  concept  of  shock  perturbation  on  external  flows  by 
using  a  short  duration  arc  discharge,34  as  seen  in  Fig.  38.  In  their  work,  the  surface  discharge  created  a 
shock-wave  perturbation  on  one  side  of  a  cone,  which  was  used  to  generate  a  steering  moment.  Future  work 
may  explore  using  nanosecond  DBD  discharges  for  an  internal  flow-path  to  achieve  similar  shock  control. 


DBD 


Figure  37.  Schematic  illustrating  a  shock  train  pertur¬ 
bation  due  to  nanosecond  DBD  actuator  locating  at  the 
leading  edge  of  a  supersonic  inlet. 


Figure  38.  Evolution  of  arc-discharge  generated  disturbance  on  the  surface  of  a  supersonic  cone,  taken  from 
Ref.  34. 

In  addition  to  exploration  of  the  controlling  shock  boundary-layer  interactions  using  nanosecond  DBD 
actuators  for  internal  flow  paths,  it  is  also  necessary  to  improve  the  fidelity  of  the  phenomenological  energy 
deposition  model.  To  this  end,  future  work  will  incorporate  the  recent  work  by  Poggie  et  al.,16  which  more 
accurately  models  the  energy  transfer  processes  in  a  pulsed  surface  dielectric  barrier  discharge. 
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